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ABSTRACT 

In two previous papers we formulated and solved, for a fiducial set of free parameters, the problem 
of the formation and evolution of a magnetically supercritical core inside a magnetically subcritical 
parent cloud. The evolution was followed into the opaque phase that resulted in the formation of 
a hydrostatic protostellar core. In this paper we present a parameter study to assess the sensitivity 
of the results (1) to the density at which the equation of state becomes adiabatic; (2) to the initial 
mass-to-flux ratio of the parent cloud; and (3) to ionization by radioactive decay of different nuclei 
("^^K and ^^Al) at high densities (rin > 10^^ cm~'^). We find that (1) the results depend only slightly 
on the density at which the onset of adiabaticity occurs; (2) memory of the initial mass-to-flux ratio 
is completely lost at late times; and (3) the precise source of radioactive ionization alters the degree 
of attachment of the electrons to the field lines (at high densities), and the relative importance of 
ambipolar diffusion and Ohmic dissipation in reducing the magnetic flux of the protostar. The value 
of the magnetic field at the end of the runs is insensitive to the values of the free parameters and 
in excellent agreement with meteoritic measurements of the protosolar nebula magnetic field. The 
magnetic flux problem of star formation is resolved for at least strongly magnetic newborn stars. A 
complete detachment of the magnetic field from the matter is unlikely. The formation of a "magnetic 
wall" (with an associated magnetic shock) is independent of the assumed equation of state, although 
the process is enhanced and accelerated by the formation of a central hydrostatic core. 

Subject headings: ISM: clouds - ISM: dust ~ magnetic fields - MHD - stars: formation - shock waves 



1. INTRODUCTION 

We previously formulated the problem of the 
ambipolar-diffusion-initiated core formation and evolu- 
tion in self-gravitating, magnetically supported model 
clouds (Tassis & Mouschovias 2006a, hereafter Paper I), 
and presented the solution for a typical model cloud (Tas- 
sis & Mouschovias 2006b, hereafter Paper II). In this pa- 
per, we quantify the dependence of the results on the 
values of the relevant free parameters of the problem. 

The system of MHD equations governing the evolu- 
tion of the model cloud, which we presented in Paper 
I, contains four free parameters: the dimensionless ini- 
tial central mass-to-flux ratio /id.co of the reference state 
relative to its critical value for collapse; the character- 
istic lengthscale Zj-cf of the initial column density distri- 
bution relative to the thermal critical lengthscale; the 
ratio of the external thermal pressure and "gravitational 
pressure" in the central flux tube of the initial reference 
state, Poxt; and the density at which the transition from 
an isothermal to an adiabatic equation of state occurs, 
riopq. Additional parameters involved in the chemical 
model are well constrained observationally and experi- 
mentally (see Appendix in Desch & Mouschovias 2001 
and references therein), and values from these studies 
are used in the calculations. 

Of the parameters referred to above, the value of Poxt 
has no significant effect on the formation and evolution 
of the supercritical core in a self- gravitating model cloud. 
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The value of Ziof is a measure of the strength of the initial 
thermal-pressure forces relative to gravitational forces in- 
side the model clou d. If^ref > 3 then th e cloud is ther- 
mally supercritical. iBasu fc Mouschovias I ()1995D showed 
that, as long as the cloud is significantly thermally su- 
percritical, the value of ^ref plays no role other than to 
determine the size and hence the mass of the initial cloud 
(which they found to be proportional to /^^f). Hence, the 
effect of varying Pcxt and /rcf is not examined further 
in this parameter study. The values used for these pa- 
rameters in all the models presented in this paper are 
Poxt = 0.1 and Z~ef = 5.57r = 4.07 x (3\/2). 

The free parameters varied are the density at which the 
equation of state becomes adiabatic, Uopq, and the initial 
mass-to-flux ratio of the parent cloud /id in units of the 
critical value. We also examine how the results change 
when, at high densities (when natural radioactivity dom- 
inates cosmic rays as the primary ionization mechanism) , 
radioactive ^^Al rather than ^'^K is the most abundant 
decaying element. 

Table [T] shows the relevant parameters of the models 
studied in this paper. Model / is our fiducial run that 
was presented in Paper II. In § [2] we examine how the 
results of the calculations change as the density riopq, at 
which the transition from isothermality to adiabaticity 
occurs, is varied. For this purpose, we compare results 
from three models: nlO, / and nl2, in which all pa- 
rameters have the same values as in the fiducial run, 
except riopq; it is equal to 10"^°, 10^^ and 10^^ cm~'^ 
in the three models, respectively. These values cover 
the expected density range in which the collapsing core 
will become optically thick (iGaustad I [l9 63: Hav ashil 
[1961 iDraine fc Leel [l98l lOssenkopf fcHennind 11991 
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TABLE 1 
Model Parameters 



Model 


riopq (cm ^) 


/Jd.cO 


dominant radionuclide 


/ 


10" 


0.25 




nlO 


10"' 


0.25 


4UK 


nl2 


1012 


0.25 


40K 


i 


always isothermal 


0.25 


40k 


m5 


lO^i 


0.50 


40 K 


m7 


10" 


0.70 


40k 


m9 


10" 


0.90 


40k 


Al 


10" 


0.90 


26 Al 



iMasunaga fc Inutsuka I Il999f ) . As part of this parame- 
ter study, we also present in § 12.11 results for model i, 
a "control run" , in which the isothermal approximation 
is retained throughout the calculation. Although clearly 
unrealistic, this run is used to assess which effects are 
caused by the change in the equation of state. The effect 
of the initial mass-to-flux ratio on the solution is dis- 
cussed in § |3l We present results for four models, /, to5, 
m7 and m9, with parameters always the same as those of 
the fiducial run, except /id,cO, which is equal to 0.25, 0.5, 
0.7 and 0.9, respectively. Finally, in § |4]we examine how 
the (late-time) results depend on the kind of radioactive 
element that dominates the ionization process at high 
densities. If the mass-to-flux ratio of the parent cloud is 
close to critical, then the evolution is much more rapid 
than in the fiducial model, and radionuclides with half- 
lives shorter than may contribute significantly to the 
ionization at high densities. For this reason, we compare 
models m9 and Al, both of which have fJ.d,cO = 0.9, but 
assume different dominant radionuclides (^°K and ^^Al, 
respectively) . 

2. DEPENDENCE ON Nqpq 

In this section we examine the effect of varying the 
density at which the assumed equation of state changes 
from isothermal to adiabatic. We compare three models: 
nlO (nopq = 10^° cm-3); / (riopq = 10" cm'^); and nl2 
{riopq = 10^^ cm~^). Model / is the fiducial run, which 
was discussed in detail in Paper II. All three models have 
an initial mass-to-flux ratio //d.co = 0.25. 

Figure ^ exhibits central quantities of the model 
clouds as functions of the central number density of neu- 
trals. In all cases, the solid line corresponds to model 
/, while the dashed and dashed-dot lines correspond to 
models nlO and nl2, respectively. 

Figure ((T^) shows the central temperature (measured 
in K) of the three models as a function of the central 
neutral density. At low densities, while the equation of 
state is still isothermal, the central temperature remains 
constant and equal to the 10 K temperature of the refer- 
ence state. In each model, at a density equal to nopq, the 
equation of state becomes adiabatic and the temperature 
begins to deviate from its initial value. The evolution of 
temperature follows the adiabatic law, T oc p'^~^. The 
break in the power law and the subsequent, less steep 
increase, is the result of the excitation of the rotational 
degrees of freedom of the neutral molecules, which cor- 
responds to a shift in the adiabatic index from 5/3 to 
7/5. In the models we implement this change at a tem- 
perature of 200 K. The models do not all reach the same 



central density when each run is stopped. This is because 
the termination criterion for each run is that the central 
temperature reaches 1000 K, at which thermal ionization 
(not accounted for in the model) is expected to become 
important. This temperature is reached at a different 
central density in each model (6 x 10^^, 6 x 10^"*, and 
6 X 10^^ cm~^ in models nlO, /, and nl2, respectively), 
depending only on the value of nopq. 

Figure ((TJd) displays the dependence of the central 
mass-to-flux ratio, normalized to its critical value, on 
the central density of neutrals. The behavior in all cases 
is very similar to that of the fiducial run (model /). The 
three curves deviate slightly from one another close to the 
density riopq of each model, while they quickly asymptote 
to an identical behavior at higher nn,c- This is because 
Bz exhibits only a mild dependence on the riopq at high 
densities. 

The evolution of the central magnetic field with den- 
sity (both z— and r— components, shown in Figs. [TJ; and 
[TJi, respectively) exhibits a relatively weak dependence 
on riopq. The most pronounced difference is caused by 
a steepening of the density dependence of Bz^c immedi- 
ately after the onset of adiabaticity. If the equation of 
state had remained isothermal, Bz,c would have turned 
over and would have asymptoted to a constant value at 
high densities (see Fig. [7h below). The behavior of Bz^c 
is similar in all three models after the onset of adiabatic- 
ity. However, in model nl2 the overall value of Bz^c is 
smaller at a given density than that in models nlO and /. 
This is so because in model nl2, adiabaticity is turned 
on at a density greater than that at which the profile 
of Bz.c begins to flatten, while in the other two models 
adiabaticity is turned on at a density smaller or about e 
qual to that at which the flattening would begin. Despite 
this fact, all runs end at a similar value of Bz^c ~ 10^/iG, 
because the final central density of each model increases 
with riopq. 

The r— component of the magnetic field at the center 
of the contracting core _Br,c reaches a maximum at a den- 
sity nn,c ~ 10^^ cm~^ and decreases beyond that. This 
is because at that density the magnetic field begins to 
decouple from the matter and the field-line deformation 
can no longer be sustained. 

Figure ([2]) shows the radial dependence of the density 
of the neutrals, rin, of its logarithmic derivative, s = 
dlnnn/dr, and of the temperature, T, at different times. 
The top, middle and bottom rows refer to models rilO, 
/, and nl2, respectively. Each curve corresponds to a 
central density enhancement of one order of magnitude 
with respect to the previous curve, except for the last 
curve, which always corresponds to a central temperature 
of 1000 K. The "star" on a curve marks the location of 
the supercritical core boundary. 

The rightmost column of Figure [2] shows profiles of the 
temperature. Only the five last curves in each case show 
deviation from the reference-state temperature. The ra- 
dius at which adiabaticity sets in does not change ap- 
preciably with time. This is so because adiabaticity is 
implemented using a number-density criterion, and does 
not change in time (see left column plots). This radius is 
somewhat closer to the center the greater riopq is, because 
the density decreases monotonically with r. The qualita- 
tive behavior of the number density of the neutrals (left 
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Fig. 1. — Central quantities as functions of the central number density of neutrals (in cm for models / (solid line, riopq = 10^^ cm '^), 



nlO (dashed line, no 



10 cm '^), and nl2 (dashed-dot line, rio 



10^ 



(a) Central temperature (in K); (b) central mass-to- 



flux ratio normalized to the critical value; (c) 2— component of the magnetic flcld (in ^G) at the center; (d) r— component of the magnetic 
fleld (in ^G) at the center. 



column) is the same as that discussed in Paper II (Fig. 
3c) for the fiducial run, with the power-law density pro- 
file steepening after the onset of adiabaticity. The final 
central density achieved by the time the temperature of 
1000 K has been reached increases with the value of riopq, 
and is therefore greatest in model nl2. 

In Paper II (Fig. 3d) we discussed the qualitative be- 
havior of the logarithmic slope of the number density 
s and pointed out that the two innermost dips are of 
thermal and magnetic origin, with the inner dip corre- 
sponding to a balance between the thermal-pressure and 
gravitational forces, and the outer dip corresponding to a 
balance between magnetic and gravitational forces. The 
overall magnitude of the thermal dip is essentially un- 
changed for different values of riopq. However, its loca- 
tion moves inward for greater Uopq as all of the post- 
adiabaticity features are translated toward smaller radii. 
By contrast, the amplitude of the magnetic dip does 
change with riopq! it is greater for smaller riopq- 



Figure [3] shows radial profiles of the z— component 
of the magnetic field (left column) and of its spatial 
derivative, h — dlnBz/dlnr (right column). All mod- 
els exhibit similar behavior. There is a fiat inner region 
followed by a local maximum at the boundary of the 
hydrostatic core. Just outside, the magnetic field de- 
creases rapidly and the profile becomes steeper because 
the transition to adiabaticity occurs at smaller densi- 
ties. In the isothermal regime there is a more grad- 
ual decrease with radius (mean slope ~ —0.7) out to 
the boundary of the supercritical core. The profile of 
the magnetic field at the cloud envelope is not affected 
by the evolution of the supercritical core. Finally, the 
break in b moves from r/Ro « 2 x 10~^ in model nlO to 
r/Ro 8 X 10^^ in model nl2, and the maximum value 
oi Bz from r/Ro w 10^'"' in model nlO to v/Rq w 5x 10^^ 
in model nl2. 

Figure [3] exhibits radial profiles of the ratio of thermal- 
pressure and gravitational forces (left column) and mag- 
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Fig. 2. — Radial profiles, at different times, for models nlO (top row, rtopq = lO'^'^ cm ^), / (middle row, rio 



lO" cm-3), and nl2 



(bottom row, rio 



10' 



The times are chosen such that neighboring curves differ by a factor of 10 in central neutral density. 



Left column: central number density of neutrals normalized to its initial value in the initial equilibrium state (rin^cO 
column: spatial derivative of the density s = dlnnn/cJln r; right column: temperature (in K). 



: 2711 cm-^); middle 



netic and gravitational force (right column) at different 
times. Each curve corresponds to a central density en- 
hancement of one order of magnitude with respect to 
the previous curve, except for the last curve, which al- 
ways corresponds to a central temperature of lOOOK. The 
thermal-pressure force becomes comparable to the grav- 
itational force in the opaque core only after the onset of 
adiabaticity. The location of the thermal shock coincides 
with the hydrostatic core boundary. Because the density 
decreases monotonically with radius and the extent of 
the opaque region is limited by the condition n > nopq, 
the size of the hydrostatic core decreases with increasing 

^opq- 

The ratio of magnetic and gravitational forces ex- 
hibits the local maximum discussed in Paper II, which 
is the cause of the magnetic shock, and the correspond- 
ing magnetic feature in the density profile. As n^pq in- 
creases, the location of this maximum is displaced to- 
ward smaller radii, and the ratio Fm/Fq in the hydro- 



static core reaches smaller values. This is so because, as 
riopq increases, the transition from isothermality to adia- 
baticity takes place at higher densities and the termina- 
tion condition for each run (T = 1000 K) is achieved at 
greater central densities. 

The small-amplitude oscillations which appear in the 
ratio of magnetic and gravitational forces inside the hy- 
drostatic core are numerical noise due teh unavoidable 
inaccuracies involved in calculating spatial derivatives 
(needed to derive the magnetic force) inside the hydro- 
static core, where the magnetic field is very nearly spa- 
tially uniform. However, because in this region the mag- 
netic field is dynamically unimportant, the impact of 
such noise on our results is insignificant. The discon- 
nected line in panel (f) is due to logarithmic plotting (in 
the region where the line has not been plotted, the value 
of the magnetic force has become numerically zero). 

Figure [5] shows the spatial profiles of the radial veloci- 
ties of the neutrals (left column) and the electrons (right 
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Fig. 3. — Radial profiles of the 2:— component of the magnetic field, normalized to the field strength of the reference state (left column), 
and of its spatial derivative, b = dlniJz/dlnr (right column) at different times, as in Figure (|2j. Top row: model nlO; middle row: model 
/; bottom row: model nl2. 
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Fig. 4. — Radial profiles of the ratio of thermal-pressure and gravitational forces (left column), and magnetic and gravitational forces 
(right column) at different times, chosen as in Figure (|2j. Top row: model nlO; middle row: model /; bottom row: model nl2. The "star" 
marks the location of the supercritical core boundary. 
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column) in units of the sound speed (in the neutrals) 
of the reference state, for the three models nlO, / and 
nl2 (top to bottom row). All models exhibit the same 
qualitative behavior: There are always two features in 
the velocity of the neutrals, the "thermal" shock and the 
"magnetic" shock at somewhat greater radii. The ve- 
locity of the electrons (which are, except at very small 
radii in model nl2, well attached to the magnetic field) 
exhibits only the magnetic 

shock. Its effect is much more pronounced than in the 
case of the neutrals, and the electrons decelerate almost 
to a halt, at the boundary of the hydrostatic core. 

The positions of the magnetic and the thermal shocks 
move to smaller radii with increasing riopq, because the 
boundary of the hydrostatic core (which moves inward 
with increasing nopq) determines the position of the ther- 
mal shock and, by extension, the position of the magnetic 
feature, which is caused by field lines piling up just out- 
side the hydrostatic core. The maximum infall velocity of 
both neutrals and electrons increases somewhat as /lopq 
increases, because the delay of the onset of adiabaticity 
implies accelerated infall for longer times. In addition, 
the reacceleration of the neutrals between the magnetic 
and the thermal shocks is greater with increasing riopq. 
This is because, as seen in FigurelU in the case of smaller 
riopq the thermal-pressure force as a fraction of the grav- 
itational force increases faster with decreasing radius be- 
tween the magnetic shock and the hydrostatic core, and 
hence the reacceleration is milder. 

The left column of Figure [6] provides crucial informa- 
tion concerning the importance of the magnetic field dur- 
ing the late stages of star formation and whether a com- 
plete detachment of the field from the matter can occur 
before thermal ionization reattaches matter to the field 
lines at T « 1000 K. At the late stages of contraction 
of the supercritical core, which are for the first time fol- 
lowed in this series of papers, ions have detached from 
the magnetic field lines and follow the motion of the neu- 
trals. Hence, the electrons are the only particles that may 
remain attached to the field lines until the onset of ther- 
mal reionization of the core. The degree of attachment 
of electrons to the field lines is measured, as we have dis- 
cussed in Papers I and II, by the attachment parameter 
Ac. The value Ac « 1 indicates that the electrons are 
well coupled to the field lines, while Ac « indicates 
that the electrons have detached and follow the motion 
of the neutrals. 

As seen in the left-column plots of Figure [6l in the 
model with Uopq = 10^°cm~'^, the electrons remain well 
attached to the magnetic field throughout the run. In 
the case riopq = 10^^ cm~'^, the electrons begin to detach 
before the central temperature reaches 1000 K, at which 
the run is stopped. However, their degree of attachment 
is still high. Finally, in the model with riopq = 10^^ cm~'^, 
the electrons in the core have detached almost completely 
before the end of the run (Ac ~ 0.1). 

Using the indirect attachment parameter of the elec- 
trons to trace the densities at which the electron de- 
tachment begins and comes to completion, we find that 
electrons begin to detach at about 10^** cm^^ and detach 
completely by a density « 10^^ cm~'^, regardless of the 
value of nopq. Hence, we arrive at the conclusion that 
the detachment of electrons ( and consequently the com- 
plete decoupling of matter) from the magnetic field lines 



takes place at a number density « 10^^ cm"^, provided that 
this density can be reached before thermal ionization be- 
comes important. In order to confidently assert whether 
the required density is achieved before thermal ioniza- 
tion sets in, and hence whether this complete decoupling 
of matter from the magnetic field occurs in nature, the 
temperature history of the collapsing core would need to 
be followed in detail, accounting for radiative transfer. 
However, given the results of our parameter study, the 
detachment of electrons before the onset of thermal ion- 
ization is unlikely. This is so because the actual temper- 
ature will increase more by the time the central density 
becomes 10^^ cm'^ than it does in the isothermal approx- 
imation up to 10^^ cm~^, beyond which adiabaticity was 
assumed. 

A second important issue during the very late stages 
of collapse is whether and when Ohmic dissipation over- 
takes ambipolar diffusion as the dominant flux reduc- 
tion mechanism. This issue is explored by the right- 
column plots of Figure m which show radial profiles of the 
ratio of the ambipolar-diffusion and Ohmic-dissipation 
timescales. The ambipolar-diffusion timescale is smaller 
at all radii and times, except in the very late stages of 
the contraction, at densities > lO^'^cm^'^. For densities 
~ 10^® cm^'^, the Ohmic-dissipation timescale becomes 
w 100 times smaller than the corresponding ambipolar- 
diffusion timescale. Thus the density at which Ohmic dis- 
sipation becomes more effective than ambipolar diffusion 
is not sensitive to the thermal history of the protostel- 
lar core, but the size of the region in which this happens 
is: The greater the density at which the protostellar core 
becomes opaque, the smaller the region in which Ohmic 
dissipation dominates. 

2.1. The Isothermal Control Run 

In this section we present results of a "control run" 
for which the assumption of isothcrmality is retained 
throughout the evolution of the contracting model cloud. 
These results complement the parameter study with re- 
gard to Tiopq, by allowing one to identify the features 
that can be attributed to the transition between equa- 
tions of state. All parameters other than the equation of 
state are identical with those of the runs discussed in the 
previous section. 

Figure ^ shows the evolution of central quantities for 
this model. The dependence of the central density of 
neutral particles on time is plotted in Figure ([7^). This 
curve is essentially identical to the corresponding one of 
the fiducial run discussed in Paper II, since the slow, 
ambipolar-diffusion-controUed early stage of the forma- 
tion of the supercritical core occurs isothermally in both 
models. 

The evolution of the central mass-to-flux ratio and of 
the r-component of the magnetic field at the center of 
the cloud, shown in Figures ([Tb) and ([111), respectively, 
remained essentially the same for different values of riopq, 
as we saw in the previous section. It is also very similar 
in the isothermal run. However, the evolution of the z- 
component of the magnetic field, shown in Figure ^) , 
shows a qualitatively different behavior at late times. In 
the previous section it was shown that, after adiabatic- 
ity sets in, the dependence of B^^c on density steepens, 
and Bz^c keeps increasing until the end of the run. In 
the isothermal model, this steepening never takes place. 
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Fig. 5. — Spatial profiles of the radial velocities of the neutrals (left column) and the electrons (right column), normalized to the sound 
speed of the reference state, at different times chosen as in Figure Top row. model nlO; middle row. model /; bottom row. model nl2. 
The "star" marks the location of the supercritical core boundary. 
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Fig. 6. — Radial profiles of the attachment parameter of the electrons (left column) and the ratio of ambipolar-difi^usion and Ohmic- 
dissipation timescales (right column) at different times, chosen as in Figure ((2J. Top row. model nlO; middle row. model /; bottom row. 
model nl2. The "star" marks the location of the supercritical core boundary. 
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Instead, the curve flattens and B^.c asymptotically ap- 
proaches 0.28 G, a value which is about 800 times the 
magnetic field of the initial reference state, and slightly 
greater than the one reached in the adiabatic runs by the 
time the temperature became 1000 K. 

Figures ([It) and ([8ji) show the neutral number density 
profile and its derivative, s = d In nn/c? In r, respectively, 
at different times chosen as in Figure ([2]). As before, each 
curve corresponds to an order of magnitude enhancement 
in the central density relative to the previous curve. In 
§ [21 (in relation to Fig. [2]) and in Paper II, we iden- 
tified two types of features inside the supercritical core 
but outside the innermost uniform-density region in the 
density profile and its derivative: an inner, "thermal" 
steepening of the profile, corresponding to a dip in the 
slope s, and a second, "magnetic" steepening, due to 
the local maximum in the z— component of the magnetic 
field and the corresponding balance between magnetic 
and gravitational forces. In the isothermal run, the in- 
ner thermal feature is completely absent, because the 
temperature remains always at 10 K and a hydrostatic 
core never forms. The magnetic feature is still present, 
but is less pronounced because the enhanced "pile-up" of 
matter and magnetic flux outside the hydrostatic core is 
now absent. 

Spatial proflles of the velocities of the neutral molecules 
and the fleld lines are shown in Figure ((8^) and ([8|d), re- 
spectively. Again, only the magnetic feature is present 
here, at r = 10~^Rq, and it is much weaker than the cor- 
responding magnetic shocks seen in the runs that allow 
adiabaticity at the highest densities. The reason for this 
difference is that the presence of the hydrostatic core 
leads to a more pronounced accumulation of flux just 
outside its boundary. 

Figures ^) and ([Sf) display spatial profiles of and 
its logarithmic derivative, 5, at different times as in Fig- 
ure Q . Although the qualitative behavior of B^ is sim- 
ilar to that of the models discussed above (moving from 
larger to smaller radii, a steepening of the power-law 
profile, followed by a local maximum and then the in- 
nermost flat region), all features are less pronounced in 
the present, isothermal case. This is especially evident 
in the case of the local maximum, because of the "pile- 
up" of field lines and matter that occurs just outside the 
hydrostatic core in the adiabatic case. 

Figures ([9^) and ([9)3) show the spatial variation of the 
ratios of thermal-pressure and gravitational forces, and 
magnetic and gravitational forces, respectively. In the 
isothermal model, the thermal-pressure force never bal- 
ances the gravitational force. However, in the innermost 
part of the collapsing supercritical core at late times, 
the thermal-pressure force overtakes the magnetic force. 
In this region, then, at late times, the thermal-pressure 
forces provide the main opposition to gravity even in this 
isothermal run. 

The spatial variation of the attachment parameter of 
the electrons, Ac, is shown in Figure ([Qt) at different 
times, as in Figure ([2|). As discussed above, the detach- 
ment of the electrons is primarily density related. The 
electrons are completely detached from the field lines by 
rin ~ 10^^ cm-3. 

The ratio of the ambipolar-diffusion and Ohmic- 
dissipation timescales is shown in Figure [9ji as function 
of radius, normalized as usual to the initial cloud radius 



Rq, at different times. Below n^^c ~ 10^** cm'^, ambipo- 
lar diffusion is much more effective than Ohmic dissipa- 
tion in increasing the mass-to-flux ratio in any given flux 
tube. It is only for densities greater than « 10^'* cm'^ 
that Ohmic dissipation becomes more important than 
ambipolar diffusion. This density is slightly greater than 
in the adiabatic cases. 

Having examined the limiting isothermal case, we con- 
clude that the results concerning the structure of the 
magnetic field, the evolution of the ratio of ambipolar- 
diffusion and Ohmic-dissipation timescales, and the den- 
sity at which the electrons detach from the magnetic 
field, are not sensitive to the assumed equation of state. 
However, the presence of a hydrostatic core in the adia- 
batic models accelerates the formation of the "magnetic 
wall" (the evolution of which eventually leads to the se- 
ries of magnetically driven shocks studied in Tassis & 
Mouschovias 2005a, b). 

3. DEPENDENCE ON THE MASS-TO-FLUX RATIO 

In order to quantify the effect of the initial mass-to- 
flux ratio of the parent cloud on the evolution of the 
supercritical core, we compare four models (/, m5, m7 
and m9). All have the same value of riopq = 10^^ cm^"^, 
but an initial central mass-to-flux ratio (in units of the 
critical value for collapse) equal to 0.25, 0.5, 0.7 and 0.9, 
respectively. In Figures pU]) and pT|) . model / is repre- 
sented by a solid line, model m5 by a triple-dot-dashed 
line, model m7 by a dashed line, and model m9 by a 
dashcd-dot line. 

Figure pUb ) shows the central mass-to-flux ratio in 
units of the critical value as a function of the central 
number density of neutrals. Although the initial mass- 
to-flux ratio is different in these models, all memory of 
the initial conditions is lost after the formation of the 
supercritical core , and the increase of the central mass- 
to-flux ratio with density is practically identical in all 
models. A very similar behavior is exhibited by the (z— 
and r— ) components of the central magnetic field, shown 
as functions of central density in Figures [TOb andfTOH. re- 
spectively. The final values of Bz (and Br) in all models 
differ by less than a factor of 2. This behavior is even 
more pronounced in the case of Ag, the attachment pa- 
rameter of the electrons, which is shown as a function 
of the central number density of neutrals in Figure llOb . 
Ac starts out at the same value in all models (Ac — 1, 
which corresponds to perfect attachment of the electrons 
to the field lines), and evolves in the same manner in all 
models. At the end of each run, by the time a central 
temperature of 1000 K is reached, its values in the four 
models are: 0.92 (model /), 0.94 (model m5), and 0.95 
(models uriJ and to9). 

The quantity that is significantly affected by the ini- 
tial mass-to-flux ratio is the time required for a particular 
central density to be reached (see Fig. [TUb). The cloud 
spends most of the time in the early phase, until a super- 
critical core forms. Thereafter the contraction becomes 
dynamic (accelerated, but much slower than free fall). 
The total evolution time, from the initial state until the 
end of the run, ranges from 5 Myr (model to9) to 17 Myr 
(model /). 

Figure (HUf ) shows the ratio of the ambipolar-diffusion 
and Ohmic-dissipation timescales at the center of the 
cloud, as a function of the central number density of 
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value, as a function of the central density enhancement; (c) z— component of the magnetic field, normalized to the field strength of the 
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of the reference state, as a function of the central density enhancement. The "star" marks the formation of the supercritical core. 



neutral particles. Although the quantitative differences 
among models before the onset of adiabaticity are appre- 
ciable, the qualitative behavior is similar in all models. 
Models with smaller initial mass-to-flux ratio exhibit sys- 
tematically smaller tad / '''od values before the onset of 
adiabaticity. After adiabaticity sets in, the evolution of 
the timescales ratio is similar in all models. Ambipo- 
lar diffusion dominates Ohmic dissipation early on, but 
the two mechanisms become equally effective at central 
densities w 10^^ cm~^ in all cases. 

Figure pT|) shows radial profiles of different quantities 
for models /, m5, ml and m9. All curves correspond 
to that instant in time when the central number density 
of neutrals reaches the value nn,c = 10^'* cm~^. On each 
curve in these plots, the "star" marks the location of the 
boundary of the supercritical core in the corresponding 
model. 

Figure (jllb .) shows the central number density of neu- 
trals as a function of radius for each of the four models. 



At the late time corresponding to these curves, the estab- 
lished density profile is practically identical in all mod- 
els, except near the boundary of the supercritical core, 
where small deviations between models can be identified. 
Just outside the boundary of the supercritical core, the 
density profile is steeper for smaller mass-to-flux ratios, 
because the envelopes are better supported by magnetic 
forces. 

Figure (fTTb ) shows the profiles of the z— component of 
the magnetic field for the four models. The field strength 
in the initial state of the parent cloud, which is preserved 
in the outermost regions of the model cloud, is different in 
each model, and is greater in models with smaller mass- 
to-flux ratio. However, within one order of magnitude in 
radius inward from the location of the supercritical core, 
the profile of is almost the same in all models. 

Figure (fTTb ) shows profiles of the radial velocity of the 
neutral particles in units of the sound speed of the ref- 
erence state. Although the behavior is qualitatively the 



12 




Fig. 8. — Spatial profiles of physical quantities of model i (isothermal run) at different times, chosen as in Figure (a) Velocity 

of neutrals normalized to the sound speed of the reference state; (b) velocity of the field lines normalized to the sound speed of the 
reference state; (c) number density of the neutrals normalized to the initial central density of the reference state; (d) slope of the density 
profile s = dlnrin/dlnr; (e) 2— component of the magnetic field normalized to the field strength of the reference state; (f) slope of Bz, 
b = din Bz/dlnr. The "star" marks the location of the supercritical core boundary. 



13 



fa 



10 



10 







10 



10 




10 °10 "lO *10 ^lO" 



r / Ro 



~R ~fi ~4 ~2 
10 °10 ^10 ^10 ^10^ 

r / Ro 



1.0 

0.8 

oj 0.6 

0.4 

0.2 
0.0 




— 8 ~fi —4 —2 
10 10 10 10 10 



r / R 







o 
o 




10 ^10 ^10 ^10 ^10° 



r / Ro 



Fig. 9. — Radial profiles of physical quantities of model i (isothermal run) at different times, chosen as in Figure (|2j. (a) Ratio of 
thermal-pressure and gravitational forces; (b) ratio of magnetic and gravitational forces; (c) attachment parameter of electrons; (d) ratio 
of ambipolar-difTusion and Ohmic-dissipation timescales (the dotted line corresponds to tad = tqd)- The "star" marks the location of the 
supercritical core boundary. 



14 



10" 

10^ 

10^ 

10 
10" 

lO'^ 
10" 

10^ 

10^ 

10 

10 
10 







a 












A^ 






10^ 



10" 



10^ 



10" 



10*= 



10^ 



c 

























10*10^10^10'°L0^^0'* 



5 10 15 20 
time (Myr) 







: d 












■ / 





1.10 " 

1.05 

1.00 - 

0.95 

0.90 _ 





10^ 




loO 


Q 

O 


-2 
10 ^ 


t- 




~\ 
Q 
< 


10"* 


h 






10"^ 




10"^ 



10* 10® lO^lO^^lO^^lO^* 



f 








f 

,# 















lo^io^io^'io'lo'^o'* 



1 0* 10^1 0^1 o'lo'lo'* 



10* 10® lO^lO'^lO^^O'* 



Fig. 10. — Evolution of central quantities for models / (solid line, /ia^cO = 0.25), m5 (triple-dot-dashed line, /id,cO = 0-5)> (dashed 
line, /id,cO = 0-7) and m9 (dashed-dot line, /ia^cO = 0.9). (a) Central mass-to-flux ratio, normalized to the critical value, as a function of 
the central number density of neutrals; (b) central number density of neutrals as a function of time; (c) and (d) z— and r— components, 
respectively, of the magnetic field, as functions of neutral density; (e) attachment parameter of the electrons; (f) ratio of ambipolar-diffusion 
and Ohmic-dissipation timescales. 




^10 ^ 10° 



10 ^10 ^10 *10 ^10° 



10 ^10 ^10 *10 ^10*^ 



10 ^10 ^10 %0 ^10° 



/ Ro 



/ Ro 



/ Ro 



Fig. 11. — Radial profiles for models / (solid line), m5 (triple-dot-dashed line), m7 (dashed line) and m9 (dashed-dot line), when the 
central density is rin.c = 10^* cm~^. (a) Central density of neutrals; (b) z— component of the magnetic field; (c) radial velocity of neutrals; 
(d) radial velocity of magnetic field lines; (e) ratio of thermal-pressure and gravitational forces; (f) ratio of magnetic and gravitational 
forces. The "star" on each curve marks the location of the supercritical core boundary. 



15 



same for all models, the evolution of the velocities and 
the values at their local minima depend on the initial 
mass-to-flux ratio. Inward of the location of the super- 
critical core boundary, the magnitude of the infall veloci- 
ties achieved at a certain radius is greater for smaller ini- 
tial mass-to-flux ratio. The reason for this trend is that 
the dynamical contraction sets in at a greater central 
density enhancement as the initial mass-to-flux ratio de- 
creases. As a result, by any given central number density, 
the core will have spent less time in the dynamical phase 
in the models with smaller mass-to-flux ratio, and the 
resulting velocities will be smaller at any given radius. 
For this reason, the maximum infall velocity achieved 
by the neutrals increases as the initial mass-to-flux ratio 
increases, and is greatest in model m9. 

Inward of the magnetic shock, the neutrals are reac- 
celerated by gravity before the thermal shock, behind 
which they almost come to rest in the hydrostatic core. 
The reacceleration is smaller for smaller initial mass-to- 
flux ratios because, as we have seen, the density proflle 
is somewhat steeper at radii smaller than the radius at 
which adiabaticity sets in; hence, the thermal-pressure 
force is greater, presenting a greater opposition to grav- 
ity- 

This interpretation is verifled by Figure (|llb ). which 
shows proflles of the ratio of thermal-pressure and gravi- 
tational forces in each model. The thermal-pressure force 
is indeed a greater fraction of the gravitational force in 
the reacceleration region (inward of the magnetic shock) 
in the models with smaller mass-to-flux ratios. The mag- 
netic force, on the other hand, decreases rapidly with 
decreasing radius inward of the plateau of the magnetic 
force, where it balances gravity. In addition, the mag- 
netic force as a fraction of the gravitational force has 
almost identical behavior in all models at a given value 
of the density. This can be seen clearly in Figure pTf). 
which plots the ratio of the magnetic and gravitational 
force in each model. 

The profiles of the velocity of the field lines, Vf, in the 
four models are compared in Figure (jlld ). The velocity 
of the field lines exhibits only the magnetic- wall feature, 
because the deceleration at the magnetic shock brings 
them, as well as the attached species (the electrons), al- 
most to rest. Hence the presence of the hydrostatic core 
boundary does not affect the motion of the field lines 
significantly. The maximum infall velocity exhibited by 
the field lines is greater in models with greater initial 
mass-to-fiux ratios, because dynamical contraction oper- 
ates for a longer time. 

4. 26 AL AS AN ALTERNATIVE RADIONUCLIDE 

All models discussed above assumed that, at very high 
densities, when the central part of the collapsing core be- 
comes self-shielded from cosmic-ray ionization, the ion- 
ization rate is dominated by radioactivity (half-life 
= 1.25 Gyr). If, however, the initial mass-to-fiux ra- 
tio of the parent cloud is relatively close to critical, then 
the evolution is rapid and radionuclides with smaller life- 
times (such as ^^Al, with a half-life = 0.716 Myr) can 
dominate the ionization at high densities. Similarly, ^^Al 
can become important if the core happens to get enriched 
because of a nearby Supernova explosion. In this sec- 
tion we compare results from models m9 and Al, both of 
which have an initial central mass-to-fiux ratio equal to 



0.9 times the critical value and riopq — 10^^ cm"'^. Their 
difference lies in the radionuclide assumed to dominate 
high-density ionization, which is '^^K in model m9 and 
"^^Al in model Al. 

Figure p2|) shows chemical abundances of different 
species as fractions of the number density of neutral hy- 
drogen molecules (left column) and the ionization rate 
as a function of number density (right column) for mod- 
els m9 (top row) and Al (bottom row). In the case of 
model Al, the plateau reached by the ionization rate at 
high densities occurs at considerably greater values (more 
than three orders of magnitude) than in model m9. As 
a result, the degree of ionization and the chemical abun- 
dances of all charged species are significantly greater at 
high densities, when the ionization rate is dominated by 
radioactive decays. 

The greater degree of ionization at high densities af- 
fects the degree of attachment of the electrons on the field 
lines. Figure (I13|) shows radial profiles of the indirect 
(Ac, left column) and direct (wcTcn, right column) at- 
tachment parameters of the electrons, at different values 
of the central number density of the neutrals. The top- 
row plots correspond to model m9 and the bottom-row 
plots correspond to model Al. As we discussed in § [31 the 
attachment parameter in model m9 remains very close to 
1 throughout the run, and hence the electrons are always 
well coupled to the magnetic field. However, in model Al 
at high densities, the attachment parameter of the elec- 
trons falls below 0.8, and hence the electrons have started 
to detach by the end of the run. This effect is clearly not 
of magnetic origin, since the magnetic field (and the di- 
rect attachment parameter of the electrons at a given 
density) is essentially the same in the two models. In- 
stead, the effect is electrostatic: the increased abundance 
of positively-charged species (positively-charged grains 
and ions), all of which have completely detached from 
the field lines and are following the motion of the neu- 
trals by the densities of interest, exerts attractive elec- 
trostatic forces on the electrons, which thereby begin to 
detach earlier than they otherwise would. 

A second important effect of the higher degree of ion- 
ization in model Al is the difference in the dependence 
on number density of the ratio of the ambipolar-diffusion 
and Ohmic-dissipation timescales, which can be seen in 
Figure (O. The increased degree of ionization and the 
better attachment of electrons to field lines allowed am- 
bipolar diffusion to operate for higher densities and hence 
its timescale to remain smaller than that of Ohmic dis- 
sipation throughout the run. 

5. CONCLUSIONS AND DISCUSSION 

In Paper II, we presented results for the ambipolar- 
difussion-initiated formation and contraction of a pro- 
tostellar fragment inside a self-gravitating, magnetically 
supported model molecular cloud, with a fiducial set of 
values for the free parameters. We followed the evolu- 
tion of the fragment into the opaque regime, through the 
formation of a hydrostatic protostellar core and the later 
onset of thermal ionization. In this paper, we a ssessed 
the sensitivity of the solution on the values of the free 
parameters; namely, the density of neutrals riopq above 
which the equation of state becomes adiabatic, the ini- 
tial central mass-to-flux ratio /id.co of the parent cloud 
(relative to its critical value for collapse), and the ra- 
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Fig. 12. — Abundances of different species (left column) and ionization rate (right column) as functions of the central neutral density 

4900 cm~^). Top row: model m9 ('"^K is the dominant radionuclide); bottom row. model Al (^^Al 



normalized to its initial value (n,i_cO 
is the dominant radionuclide). 

dionuclide that dominates the ionization at the highest 
densities (nn ^ 10"'^^ cm"^). 

We found that varying riopq did not resuh in any 
quaUtatively or quantitatively significant change in the 
structure of the magnetic field inside the hydrostatic 
core. The value of Bz ~ 0.1 G reached by the end 
of each run (when T = 1000 K) is almost the same 
in all cases and, remarkably, very close to the values 
of the p rotosolar magnetic field as measured in me- 
teorites (llevvl 119881: IStacev. Lovering. k ParrVI 119611 : 
iHerndon fc Rowe I Il974f ). The detachment of the elec- 
trons from the magnetic field lines (signifying the decou- 
pling of the magnetic field from the matter) was found 
to occur at rin « 10^^ cm~'^, provided that this den- 
sity could be reached in the central part of the hydro- 
static core before the temperature reached 1000 K. If 
the temperature of 1000 K is reached before the elec- 
trons have detached from the field lines, thermal ioniza- 
tion will recouple the field lines and the matter; hence, 
complete decoupling will never be achieved. This latter 



case was found to be the most likely outcome, because 
the electrons did not detach from the field lines unless 
adiabaticity set in only at very high neutral densities 
(> 10^^ cm~^). However, a detailed radiative transfer 
treatment is needed in order to follow the evolution of 
the core temperature more accurately and thereby arrive 
at a definitive conclusion concerning the decoupling of 
the magnetic field from the matter. 

That the structure of the magnetic field, the density at 
which the electron detachment takes place, and the rela- 
tive importance of ambipolar diffusion and Ohmic dissi- 
pation in reducing the protostcUar magnetic flux are in- 
sensitive to the assumed equation of state is also corrob- 
orated by the control isothermal run (in which the equa- 
tion of state is assumed to remain isothermal throughout 
the core and for the entire duration of the run). How- 
ever, the onset of adiabaticity and the formation of a 
central hydrostatic core did facilitate the development of 
the magnetic shock, which, in the adiabatic cases, formed 
at earlier times and was significantly stronger. 
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Fig. 13. — Spatial profiles of the indirect (left column) and direct (right column) attachment parameters of electrons, in models m9 (top 



The effect of varying the initial mass-to-flux ratio of the 
parent cloud was found to be limited to the time required 
for the formation of the magnetically supercritical core 
(which is longer for the more subcritical clouds) . Beyond 
that point, the supercritical core loses memory of the 
initial mass-to-flux ratio of the cloud, and the late-time 
results are qualitatively and quantitatively similar in all 
cases. 

Finally, changing the dominant radionuclide (from ^''K 
to ^^Al) responsible for ionization at the late stages of the 
evolution, {un > 10^^ cm"'^) 

led to the electrons beginning to detach from the mag- 
netic field lines at slightly smaller densities, because the 
more abundant, detached, infalling positive charges en- 



hanced the electrostatic attraction of the electrons. In 
addition, ambipolar diffusion remained more effective 
than Ohmic dissipation in reducing the magnetic flux 
of the protostellar core t hroughout the run. 
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